Add CRT-PMT flash classification to TPC-PMT flash-matching producer #875
Conversation
PetrilloAtWork
left a comment
PetrilloAtWork
left a comment
There was a problem hiding this comment.
Two topics for your consideration:
- recommendation to make the extraction of CRT matching optional (e.g. by specifying an empty CRT/PMT matching label): let's not force CRT on people who do not care about it (because they just don't care, because it takes too long to execute, because it was dropped from the input, because there was no CRT component in that data run, ...)
- suggestion that it be disabled by default
- suggestion to preserve the
MatchTypetype in the data product, and in most of the code
Also a couple of code style notes.
| //Find index of flash that minimizes barycenter distance along Z, or in the YZ plane | ||
| if ( fMinimizeZDistance ) { | ||
| thisFlashCenterZ = flash.ZCenter(); | ||
| thisDistance = std::abs(thisFlashCenterZ - fChargeCenterZ); | ||
| if ( thisDistance < minDistance ) { | ||
| minDistance = thisDistance; | ||
| matchIndex = m; | ||
| } | ||
| } | ||
| else { | ||
| thisFlashCenterY = flash.YCenter(); | ||
| thisFlashCenterZ = flash.ZCenter(); | ||
| thisDistance = std::hypot( (thisFlashCenterY - fChargeCenterY), (thisFlashCenterZ - fChargeCenterZ) ); | ||
| if ( thisDistance < minDistance ) { | ||
| minDistance = thisDistance; | ||
| matchIndex = m; | ||
| } |
There was a problem hiding this comment.
Extract the common code paths from the if branches:
| //Find index of flash that minimizes barycenter distance along Z, or in the YZ plane | |
| if ( fMinimizeZDistance ) { | |
| thisFlashCenterZ = flash.ZCenter(); | |
| thisDistance = std::abs(thisFlashCenterZ - fChargeCenterZ); | |
| if ( thisDistance < minDistance ) { | |
| minDistance = thisDistance; | |
| matchIndex = m; | |
| } | |
| } | |
| else { | |
| thisFlashCenterY = flash.YCenter(); | |
| thisFlashCenterZ = flash.ZCenter(); | |
| thisDistance = std::hypot( (thisFlashCenterY - fChargeCenterY), (thisFlashCenterZ - fChargeCenterZ) ); | |
| if ( thisDistance < minDistance ) { | |
| minDistance = thisDistance; | |
| matchIndex = m; | |
| } | |
| // Find index of flash that minimizes barycenter distance... | |
| double const thisFlashCenterZ = flash.ZCenter(); | |
| if ( fMinimizeZDistance ) { // ... along Z | |
| thisDistance = std::abs(thisFlashCenterZ - fChargeCenterZ); | |
| } | |
| else { // ... in the YZ plane | |
| thisDistance = std::hypot( (flash.YCenter() - fChargeCenterY), (thisFlashCenterZ - fChargeCenterZ) ); | |
| } | |
| if ( thisDistance < minDistance ) { | |
| minDistance = thisDistance; | |
| matchIndex = m; | |
| } |
For example, this helps noticing that the minimisation code is the same, and only the metric changes.
Note that the suggestion also requires removing the definition of thisFlashCenterZ and thisFlashCenterZ from above (rule: declare as close as possible to where it is used).
Even more compact, but less readable:
| //Find index of flash that minimizes barycenter distance along Z, or in the YZ plane | |
| if ( fMinimizeZDistance ) { | |
| thisFlashCenterZ = flash.ZCenter(); | |
| thisDistance = std::abs(thisFlashCenterZ - fChargeCenterZ); | |
| if ( thisDistance < minDistance ) { | |
| minDistance = thisDistance; | |
| matchIndex = m; | |
| } | |
| } | |
| else { | |
| thisFlashCenterY = flash.YCenter(); | |
| thisFlashCenterZ = flash.ZCenter(); | |
| thisDistance = std::hypot( (thisFlashCenterY - fChargeCenterY), (thisFlashCenterZ - fChargeCenterZ) ); | |
| if ( thisDistance < minDistance ) { | |
| minDistance = thisDistance; | |
| matchIndex = m; | |
| } | |
| // Find index of flash that minimizes barycenter distance... | |
| double const thisFlashCenterZ = flash.ZCenter(); | |
| double const thisDistance = fMinimizeZDistance | |
| ? std::abs(thisFlashCenterZ - fChargeCenterZ); // ... along Z | |
| : std::hypot( (flash.YCenter() - fChargeCenterY), (thisFlashCenterZ - fChargeCenterZ) ); // ... in the YZ plane | |
| if ( thisDistance < minDistance ) { | |
| minDistance = thisDistance; | |
| matchIndex = m; | |
| } |
| matchedFlashClassification = static_cast<int>(match->flashClassification); | ||
| } | ||
| else { | ||
| std::cout << "No CRT-PMT matching information for this flash." << std::endl; |
There was a problem hiding this comment.
Best practice is to send all messages to messagefacility.
This module ignores that, but it does protect from indiscriminate output:
| std::cout << "No CRT-PMT matching information for this flash." << std::endl; | |
| if (fVerbose) | |
| std::cout << "No CRT-PMT matching information for this flash." << std::endl; |
| @@ -279,6 +284,7 @@ class TPCPMTBarycenterMatchProducer : public art::EDProducer { | |||
| double fFlashCenterZ; ///< Weighted mean Z postion of hit PMTs (cm) | |||
| double fFlashWidthY; ///< Weighted standard deviation of Y postion of hit PMTs (cm) | |||
| double fFlashWidthZ; ///< Weighted standard deviation of Z postion of hit PMTs (cm) | |||
| int fFlashClassification; ///< Flash classification according to the CRT-PMT matching | |||
There was a problem hiding this comment.
Also here I suggest to preserve the original data type as long as possible:
| int fFlashClassification; ///< Flash classification according to the CRT-PMT matching | |
| sbn::crt::MatchType fFlashClassification; ///< Flash classification according to the CRT-PMT matching |
and to convert it only when writing it to a TTree (if my suggestion on SBNSoftware/sbnobj#157 is adopted, the other conversion will not be necessary).
| void updateFlashVars(art::Ptr<recob::OpFlash> flash, double firstHit); ///< Update slice-level data members with best match info | ||
| void updateMatchInfo(sbn::TPCPMTBarycenterMatch& matchInfo); ///< Update match product with slice-level data members | ||
|
|
||
| void updateFlashVars(art::Ptr<recob::OpFlash> flash, double firstHit, int matchedFlashClassification); ///< Update slice-level data members with best match info |
There was a problem hiding this comment.
I suggest to preserve the original data type as long as possible:
| void updateFlashVars(art::Ptr<recob::OpFlash> flash, double firstHit, int matchedFlashClassification); ///< Update slice-level data members with best match info | |
| void updateFlashVars(art::Ptr<recob::OpFlash> flash, double firstHit, sbn::crt::MatchType matchedFlashClassification); ///< Update slice-level data members with best match info |
| @@ -306,10 +312,12 @@ TPCPMTBarycenterMatchProducer::TPCPMTBarycenterMatchProducer(fhicl::ParameterSet | |||
| fOpFlashLabel(p.get<std::string>("OpFlashLabel")), | |||
| fPandoraLabel(p.get<std::string>("PandoraLabel")), | |||
| fTriggerLabel(p.get<std::string>("TriggerLabel")), | |||
| fCRTPMTMatchingLabel(p.get<std::string>("CRTPMTMatchingLabel")), | |||
There was a problem hiding this comment.
I recommend to put the new feature as optional, and suggest that the default behaviour disables it to make this not breaking.
A few changes are needed for that to happen.
| @@ -569,13 +586,24 @@ void TPCPMTBarycenterMatchProducer::produce(art::Event& e) | |||
| unsigned unsignedMatchIndex = matchIndex; | |||
| const art::Ptr<recob::OpFlash> flashPtr { flashHandle, unsignedMatchIndex }; | |||
|
|
|||
| //Get CRT-PMT matching classification for the matched flash | |||
| art::FindOneP<sbn::crt::CRTPMTMatching> matchPtr(flashHandle, e, fCRTPMTMatchingLabel); | |||
There was a problem hiding this comment.
art::FindOneP should not be initialised inside a loop, since it is expensive.
This one should be immediately after the retrieval of flashHandle.
To make the thing more complicate there is my other request, of making this feature optional, since a art::FindOneP variable must be always fully initialised. My typical solution to art::FindXxx is to declare it... well, optional:
auto const matchPtr = fCRTPMTMatchingLabel
? std::make_optional<art::FindOneP<sbn::crt::CRTPMTMatching>>(flashHandle, e, fCRTPMTMatchingLabel)
: std::nullopt;(matchPtr will be defined of type std::optional<art::FindOneP<sbn::crt::CRTPMTMatching>>). std::optional objects behave like "smart" pointers, so you can check if they are assigned (if (matchPtr) ...) and you access their value by dereferencing (matchPtr->at(matchIndex)).
This PR (based off of
v10_04_04) adds the flash classification from the CRT-PMT matching (via thecrtpmtlabel) to the best-matched flash in the TPC-PMT barycenter flash-matching producer. This information can be propagated to the CAFs, and enables for using the CRT-PMT matching information at the slice-level in selections. Relevant PRs will follow.It also adds an option (via
MinimizeZDistance) to choose the flash that is best-matched to the TPC charge by minimizing only the distance along the longitudinal direction. By default, the code still chooses the best-matched flash by minimizing the longitudinal-vertical distance.Associated PRs
Review
Tagging for review @francescopoppi and @PetrilloAtWork as the CRT & light gurus. Thanks!